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DESIGN AND GLOBAL ANALYSIS OF SPACECRAFT 


ATTITUDE CONTROL SYSTEMS 
George Meyer 
Ames Research Center 


SUMMARY 


A general procedure for the design and analysis of three-axis, 
large-angle attitude control systems has been developed. Properties of 
three-dimensional rotations are used to formulate a model of such systems. 

The model is general in that it is based on those properties which are common 
to all attitude control systems, rather than on special properties of partic- 
ular components. Numerical values are assigned to attitude error by means of 
error functions. These functions are used to construct asymptotically stable 
control laws. The overall (global) behavior of the system is characterized by 
the envelope of all time histories of attitude error generated by every possi- 
ble combination of initial condition, target attitude motion, and disturbance. 
A method for computing upper bounds on the response envelope is presented. 
Applications of this method indicate that it provides a useful alternative to 
Liapunov analysis for the determination of system stability, responsiveness, 
and sensitivity to disturbances, parameter variations, and target attitude 
motion . 


INTRODUCTION 


The complete design of an attitude control system, generally speaking, has 
four phases; namely, (1) specification of the task to be performed by the sys- 
tem and the selection of major system components; (2) design of the controller 
linking sensors to torquers; (3) verification of the design; and (4) construc- 
tion of the system. The present report deals primarily with phases (2) and 
(3). 


One approach to the design of a controller is provided by the theory of 
optimal control. The methods of this theory are elegant and explicit. Unfor- 
tunately, they are difficult to apply when the system is both nonlinear and 
multidimensional which is the case for large-angle, three-axis attitude con- 
trol systems being considered here. Not only is the computation of control 
laws for such systems very time consuming, but the control laws, once computed, 
are difficult to implement. These difficulties are responsible for the 
limited enthusiasm shown in the field for the routine application of optimal 
control theory to the design of control laws for complex systems. Neverthe- 
less, this theory is very useful for the analysis of system performance. 

Another approach, which is most frequently taken in practice and, also, 
the one taken here, is to pick a structure of the control law which on 



physical grounds seems most likely to be adequate, and then test the resulting 
system by means of a computer simulation. The success of this approach 
depends on two factors: (1) the familiarity of the designer with the general 

properties of the system, and (2) the validity of the inference that is made 
from the results of the simulation. The first is obvious. However, the 
second needs some elaboration, particularly since often it is either 
overlooked entirely or dismissed as insignificant. 

In a simulation one is concerned with accuracy and completeness. A 
simulation is accurate if the error between the simulated response and the 
response of the actual system to a given control situation (i.e., initial con- 
dition and forcing function) is small. A simulation is complete if the 
behavior of the system for any possible control situation can be inferred from 
the collection of cases simulated. If all possibilities are accurately simu- 
lated then the simulation is, of course, complete. However, in most cases 
occurring in practice the number of possible control situations is much larger 
than is practical to simulate. In such cases, it is necessary to justify the 
inference that is made in going from the limited collection of tests to the 
overall (global) behavior of the system. Without this, one cannot be certain 
that all possible control situations resulting in system failure have been 
included in the test sample. Thus, the ad hoc design procedure is realistic 
only if it is followed by a global analysis in which the behavior of the sys- 
tem, or at least an upper bound on the worst case, for all possible control 
situations is determined analytically. This means that in the selection of 
candidate control laws the ease with which the subsequent global analysis can 
be performed must be kept in mind. A control law that is easy to implement 
and that seems, on physical grounds, to be adequate may, in certain cases, be 
rejected if it leads to system complexity beyond the reach of the available 
techniques of global analysis. 

Attitude control systems vary considerably in internal structure. Thus, 
for example, torque may be generated by means of reaction wheels, control 
moment gyros, or reaction jets. Similarly, attitude may be measured by means 
of star trackers, sun and magnetic field sensors, or inertial gyros. Angular 
velocity may be measured directly, as with rate gyros, or it may be computed 
from attitude data. This variety has resulted in a corresponding variety of 
control laws, each using special properties of particular components (see, 
e.g., refs. 1 and 2). In most cases these designs are based on sound physical 
insight, but they result in systems that are difficult to analyze, and the 
analysis is not carried out. It is, therefore, desirable to approach the 
design problem from a general point of view by taking advantage of the basic 
similarity of all attitude control systems; namely, that (1) their primary 
control objective is to maintain the spacecraft as close to the desired atti- 
tude as is necessary for successful operation, and (2) their basic nonlinear- 
ity is that of three-dimensional rotations. This approach was taken in 
references 3 to 8, and the present report may be considered to be a 
generalization of that work. 

A general model of attitude control systems is formulated from the 
properties of three-dimensional rotations. Several representations of the 
distance between the spacecraft and target attitudes are given in terms of 
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attitude error functions. These functions are used, first, to generate 
control laws for which the system is asymptotically stable, and then they are 
used to characterize the overall performance of the system by means of its 
response envelope. This envelope is defined at each instant in the control 
interval as the maximum of the attitude error generated by every combination 
of admissible initial condition, target attitude motion, and disturbance. 

Worst case performance of the system may be estimated by means of upper bounds 
on the response envelope. Methods for computing such upper bounds are pre- 
sented, Finally, the design and analysis are illustrated with examples. 


SOME PROPERTIES OF THREE-DIMENSIONAL ROTATIONS 


The equations of motion of a spacecraft attitude control system are 
determined to a large extent by the properties of three-dimensional rotations. 
Some of these properties are summarized in the present section. 


Representation of Attitude - Attitude Error Equation 

The attitude of a rigid body may be defined relative to a given reference 
by means of a pair of right-hand orthonormal triplets of vectors with a common 
origin 0. Let the triplet s = ^s2> fixed in the given refer- 

ence space, and let the triplet a = (e^^, e^2> ^a3^ fixed in the body. In 
the usual case, the given reference will be inertial space, and the common 
origin of the triplets will coincide with the body center of mass. Let the 
transformation which maps s into a be denoted by ^as > that 

Oai ~ ^as^si^ ^ ~ ^ (^3 


and let be represented in the s-basis by the 3x3 matrix A^g. This 

matrix will henceforth be interpreted as spacecraft attitude with respect to 
inertial space. The elements a^^j of A ^^3 are the direction cosines between 

a and s, that is, aij = eai * gives the 

coordinates of e^i with respect to s. Since both s and a are right hand 
and orthonormal, A^g is a rotation matrix; that is, detCAg^g) = 1 and 


where C denotes the transpose of ( ) , and I is the identity matrix. 

Let X be an arbitrary vector, and let its coordinates in the s and a 
basis be Xg and Xg, respectively; then 


X 


a 


= A 


as^s 


( 3 ) 
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Let y be another vector 
be the cross product of 
of z are given by Zg = 
defined for any u in 


. Then the dot product x • y = x^yg = x^y^. Let z 
X and y; that is, z = x x y. Then the s coordinates 
-S(xs)ys, where the skew-symmetric matrix S is 
by 


/“ 

U 3 

-U2\ 


S(u) = I-U3 

0 


(4) 

\u2 

-Ui 

0 / 



Similarly, z^ = -S(Xa)ya- Hence, A^^SCx^DA^g = S(xg), and, in general, for 
any x in and any rotation matrix A, 


A^S(Ax)A = S(x) 


(5) 


In addition, the matrix S has the following useful properties. For any u 
in E^, 


S^(u) = -u^ul + uu^ 


( 6 ) 


and for any u and v in E^, 


S(u)v = -S(v)u 


(7) 


and 


S(u + v) = S(u) + S(v) (8) 

These relations play an important role in what follows. As an immediate 
application, consider the problem of measuring A^. Suppose that two iner- 
tially fixed directions are given by two unit vectors xl and x^. If the 
inertial coordinates Xg, x^ of these vectors are known, and if their body 
coordinates x^, x^ can be measured, then the attitude matrix A^g can be 
coii 5 )uted as follows. Form the matrix Bg in which the first two columns are 
x^ and x^ and the third is the cross product x| = -S(x^)x^. Similarly, form 
the matrix from measurements. Then it follows from equation (3) that 

Aas = ^aBg^- The inverse of Bg exists whenever the two inertial directions 
are not colinear. In an actual mechanization, the two inertial directions may 
be defined by stars, in which case x| are given by star tables, and xi may 
be computed from star tracker gimbal angles (see example 7 on p. 18). 

The target (desired) attitude may be defined by a third right-hand 
orthonormal triplet of vectors, say cl = e^^, e^^) whose origin is 

common with that of s and a. The transformation mapping s into d will be 
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denoted by Ajg . Its representation in s will be denoted by the matrix Ajs 
and identified as the target attitude. The elements of Ajg are the 
direction cosines * Cg j . 

Consider, now, a way in which attitude error may be defined. When the 
actual attitude of the spacecraft is the desired attitude, A^^^ = A^ 3 . This 
condition may be expressed either as Ag ^3 - Aj^ = or as ^as^ds ” which 
by orthogonality of is equivalent to A^^^A^^ = I, Since the ease with 

which a given set of problems can be solved depends, strongly on the under- 
lying mathematical structure, it is desirable to introduce as much structure 
at the outset as possible. The set of three-dimensional orthogonal matrices 
has a well-developed structure. In order to have this structure at hand in 
what follows, system attitude error will be defined by the orthogonal matrix 


R = 


A A^ 
^as^ds 


(9) 


It may be noted that R represents the rotation from d to a. The matrices 
Ads, Ag^^, and R will be referred to as the system input, output, and error, 
respectively. 


Angular Velocity - The Kinematic Equation 

Suppose that the spacecraft a is rotating with angular velocity oo 
relative to the inertial space s. Let x define a point P fixed in the 
body. Then, x^ = A^s^s> ^a ^ ^s Hence, the following 

chain of equations is true. 


0 = X. 


^as ^s 


AasX, 


= [A; 


as^L - AasS(A^u)JA 


as a-^ as- 




Since this chain is true for any fixed point P, it follows that. 


^as “ ^ ^^a^ ^as 


(10) 


This is the kinematic equation of spacecraft attitude. The corresponding 
equation for the target is 


^ds = (11) 

The attitude error is defined by (9) . Hence, R = A^^A^^ + A^^A^^, and from 
equations (10) and (11) it follows that 


R = S(u)a)R - RS(w^) 


( 12 ) 
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This is the kinematic equation of attitude error. Note that equation (12) is 
a well-behaved differential equation. It is defined for all attitude errors; 
it is without singularities; and, for given time histories of spacecraft and 
target angular velocities, it is linear. According to equation (5), 

RS(o)(i) = S(Ru)jj)R; hence, equation (12) is equivalent to the following equation 

which is useful for certain purposes. 


R = S(o)a - (13) 

The argument of S in (13), - Rwj = <*), can be interpreted as the a 

coordinates of angular velocity error. 

Parameterization of Attitude Error Matrix 

The elements of R are not independent; they are connected by the 
orthogonality condition RR^ = I. Thus, R may be parameterized with fewer 
than nine parameters. The minimum number is three (i.e,, Euler angles). How- 
ever, the minimum number is not always convenient. A particularly useful 
parameterization is the following. 

Suppose that a, originally coincident with d, is rotated with a 
constant angular velocity a relative to d. Then R(0) = I, and 
R(t) = S(a)R(t), where a is constant. This is a linear differential equa- 
tion with constant coefficients whose solution is R(t) = exp[S(a)t]. Let 
(j) = II a 11 1 and c = a/ II cell, so that lie II = 1. Then 

R = exp[(|)S(c) ] (14) 

The scalar 4> will be interpreted as the magnitude of attitude error, and the 
unit vector c will be interpreted as the direction of attitude error. 
Expanding the exponential in equation (14), and noting from equation (6) that 
S^(c) = -S(c), and so on for higher powers of S, and collecting terms yields 
the following equivalent expression of equation (14) : 


R = I + sin c|)S(c) + (1 - cos c())S^(c) 


(15) 


Since S(c)c = 0, c is an eigenvector of R; that is, Rc = c. The 
corresponding eigenvalue is 1 . It is the consequence of Euler's theorem on 
rotations that any attitude can be reached from I by a constant angular 
velocity. Hence, any R can be parameterized by ((fi, c) as in equation (14) 
and its equivalent equation (15) . 

The parameters (cfi, c) can be computed from the elements rj^j of R 
expressed by equation (15). Thus, since trace S(c) = 0, and trace S^(c) =-2, 
trace(R) = 3 - 2(1 - cos c()) ; hence. 
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( 16 ) 


(f) = arcos I i [trace (R) - 1] 
[0,tt] 


where [0,n] is the closed interval from 0 to ir. On the other hand, since 
both 1 and S^(c) are symmetric, sin <p S(c) = (1/2) (R - R^) . Hence, for 
0 < <}) < IT, 


J 

/r2 3 

- T32\ 


^ ~ 2 sin 1 

^31 

- ri3 j 

(17) 

\ 

Vri2 

- ^21/ 


The cases 4> = 0 ,it are singular. When 

4> = 

0, c is arbitrary. 

When 4> = TT, 


equation (15) shows that the components of c are the solutions to the 
following equations. 



1 

CiCj - y 


(17a) 


The error angle function (p defined by equation (16) has several useful 
properties, which are listed below. 

(i) In the one-dimensional case (i.e., shaft positioning servo), in which 

c is constant, the usual definition of error is <t>e ~ ‘t’d “ '^a’ 'f’d 

are input and output angles, respectively. In that case <{>= Ugl if 

Ugl i. and 4) = 2tt - |())g| if tt <. (t>g ^ 2tt . 

(ii) When 4> is small, the attitude error may be represented by the 
vector 4>c. Its components are, to first order, Euler angles of R, and 4> 
is the square root of the sum of the squares of these angles . 

(iii) Let 6j be the angle between the itk vector of a and the ith 
vector of 3. Then 0i <. <t>, because rj^j^ = cos 0^ = cos <|) + (1 - cos (P)c^^ 
which implies that cos 0j^ ^ cos (J). 

(iv) Consider all paths from I to R. Each satisfies the differential 
equation R = S[o)(t)]R for some time history m and boundary conditions 
R(0) = I, R(tf) = R. It can be shown (see appendix A) that for all m. 





IU(t)ll dt 


so that 4i may be considered to be the minimum angular distance between a 
and d. This also means that for any two orthogonal (three-dimensional) 
matrices A and B, 4)(ABt) <. 4i(A) + <ti(B). In fact, the function 
<KA,B) = (|)(ABt) is a metric on the space of three-dimensional rotation matri- 
ces: (1) <})(A,B) is positive; (2) <1)(A,B) = 0 if and only if A = B; 

(3) <|)(A,B) = 4i(B,A); and (4) 4i(B,A) + 4>(A,C) 2 . <|>(B,C) , 
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Consider, next, the kinematic equation in terms of . It is shown in 

appendix B that for any tfi, c, and m 


i = ctw (18a) 

c = J s(u))c + y cot - (w^c)c] (18b) 


Equations (18) are singular at (}> = 0 and (() = tt. To gain some insight into 
the properties of equations (18) consider the special case in which the 
angular velocity uj is a constant unit vector. Let the angle between c and 
b) be denoted by ip. Then, according to equation (18a), = cos ijj, and from 

equation (18b), (cos = (1/2) cot [(1/2) (})]sin^ p. Hence, the motion passing 
through the point ((|^q, 4 ^q) remains on the curve 



3rA 


/ \ 


% 

\\ 

\\ 


/ 


\\ 

\\ 

\ 




2tt 


\ / 

1 ^ I 
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Figure 2.— Time histories of error angle for constant 
angular velocity. 


;in(i ♦) 


sin = sin 


(l*o) 


sin 


Figure 1. Motion of system with constant angular velocity. 


Curves corresponding to three 
different initial conditions are 
sketched in figure 1. Curve is 

at the point Dj at time t = 0. As 
the body rotates about oj, it swings 
by the target. The closest approach 
is at point B where (p = 50° and c 
is perpendicular to o). Thereafter, 

4) increases until <p = 180° (point 
C) where c switches sign (point 
Di) . The body then begins to 
approach the target = 0) repeating 
the cycle. If c is colinear with 
tji) at t = 0 (point D2) , it remains 
so for all time, switching sign of c 
at 4) = 0° and p = 180° (curve F2) . 
Finally, if at t = 0, 4> = 180° and 
c is perpendicular to m (point D3) , 
then the body never approaches the 
target but remains always at the 
maximum distance away (point F3). 
Figure 2 shows the time histories of 
the error angle 4> corresponding to 
these three cases. 

The singularity at 4^ = 0 can 
be removed by multiplying c by a 
suitable function of p. For example, 
Euler parameters (n,e) are defined as 
follows . 
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(19a) 


n 


2 



e; 


2 sin 



(19b) 


These parameters are related by the equation = 4, The corresponding 

kinematic equation may be obtained from equations (18) . It is the following 
set . 


n = - y (20a) 

e = y S(o))e + y no) (20b) 


Euler parameters have the advantage that the corresponding kinematic equa- 
tions (20) is continuous at cf> = 0, and it does not involve trigonometric 
functions. However, the singularity at (f) = tt is still present. 

There are other parameterizations of three-dimensional rotations, but 
since they will not be used in this report, they will not be discussed. 


Angular Acceleration - The Dynamic Equation 

The motion of the spacecraft about the fixed point 0 may be influenced in 
two ways: By means of external torques (i.e., reaction jets); and by means of 
internal torques generated by an angular momentum exchange and storage device 
(i.e., reaction wheels, control moment gyros). In general, both influences 
may be present. Let the angular momentum of the main body (about the fixed_ 
point 0) and that stored in the^exchange device by denoted by the vectors h 
and h^, respectively, and let T be the total_extemal torque acting on the 
system. The total angular momentum h = h^ + h^, and the time derivative of 
its inertial coordinates is, according to Newton's law, 

h = Tg = T (21) 

s s as a 

Let the a-coordinates of the moment of inertia of the main body be denoted by 
the matrix J^* Then in a. 


But ha = ^as^s* according to equation (6) A^g = S(o)a)Ag^g. Hence, 
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Therefore, 


“a = + Ja'S(o>a)Aashs ' j;'ja“a ^22) 

The following interpretations will be given to the terms on the right hand 
side of (22) . The first term represents the effect of angular momentum 
exchange rate. The second term represents the effect of external torque. The 
third represents gyroscopic coupling. The fourth is due to time variation of 
body moment of inertia. Equations (21) and (22) together constitute the sys- 
tem dynamic equation. Its form is useful for control applications because the 
variables (~h^)‘, T^, and appear explicitly, and they are the ones usually 
available for control. The following special cases appear frequently in 
practice. 

Case 1 - external torque- Suppose that the angular momentum exchange 
and storage device is inactive so that h| = 0, for all t 0, and that the 
moment of inertia is constant. Then h 3 ^ = A^^h^ = dynamic 

equation is the Euler's equation of motion. 


The control variable is the external torque T^. It is typically generated by 
means of reaction jets, external magnetic field, gravity gradient, etc. 

Case 2 - reaction wheels- Suppose that the system is being controlled by 
means of an angular momentum exchange and storage device consisting of three 
reaction wheel-motor pairs placed along the body axes a. Let be the 

moment of inertia of the main body with locked wheels; let be the 

diagonal matrix whose elements are the moments of inertia of the wheels about 
their spin axes, and let Tg be the column matrix whose elements are the 
wheel motor torques. Then defining Ja “ '^a “ assuming no external 

torque and constant Ja, it follows that 

<^a = + j;'s(a>a)Aash3 (24) 

where hg is a constant. In this case the control variable is the motor 
torque T^. 

Case 3 - contvot moment gyros- Suppose that the system is being 
controlled by means of a set of control moment gyros. Let the active gimbal 
angles of all the gyros in the package be arranged in a column matrix q of 
an appropriate dimension, and let the a-coordinates of the total angular 
momentum of all gyros be denoted by h^. Assuming that the total angular 
momentum of each gyro may be adequately represented by its spin momentum, we 
may express h*' as a function of q, say, h^ = h(q) . Then (h 3 )‘= h (q)q, and 

3. 3 “ H 
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for Tji = 0 and constant the dynamic equation is given by 

“a “ l-“hq(q) ]q + SCoiaJAaghg 

where hg is a constant. Suppose that the gimbal angles are driven through a 
processor so that 

q = -F(q)Jae 

where e is the input to the processor. Then 

<^a “ "^a ^^‘^a^'^as^s (25) 


In this case e is the control variable. If h is a one to one mapping from 
the region of interest Q onto h(Q) and the processor is such that 
F(q) = hZ^Cq), then (25) acquires a particularly simple form. 

Equations (11), (13), (21), and (22) describe the system. For given 
initial condition and time histories of the control variables, the motion of 
the system is the corresponding solution of these equations. The detailed 
mathematical model is considered in the following section. 


Mathematical Model of Attitude Control Systems 

The discussion in the preceding section motivates the mathematical model 
of attitude control systems shown in table 1. The model might appear at first 
sight unnecessarily detailed; however, it is basically quite simple and useful 
for the purposes of the following discussion. The detail given is required by 
the analytical techniques to be employed. 


TABLE 1.- MATHEMATICAL MODEL OF ATTITUDE CONTROL SYSTEMS 


State space 

X = X = (R, 

^a) 

Region of operation 

0 = {xrRR't = I and 

det(R) = 1} 

Admissible initial conditions 

0o c 0 


State equation 

X = g[t,x,u(t) ,y] 


Kinematic equation 

R = S(u)a + yi)R 


Dynamic equation 

o)a = f(t,x,y,y2) 


Nominal control law 

f(t,x,y,0) 


Mode variable 

y 


Perturb ations 



Target velocity 

y^ = ni(t,x,y)ui(t) 


Disturbance 

= n2(t,x,y)u2(t) 


Admissible forcing functions 

(Ui,U2) G u 
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The underlying state space X is 12-dimensional with the state x 
denoted mnemonically by (R,o)a) . The first nine components of x are the 
elements of the error matrix R, and the remaining three coordinates are the 
body (a) coordinates of body angular velocity 0)3. The region of operation 0 
of the system is the 6-dimensional manifold defined by the condition that R 
be a rotation matrix. All possible motions of the system are inherently 
restricted to this 0 . The reason for imbedding the system in is to 

have a nonsingular kinematic equation. The set of admissible initial condi- 
tions 0Q, which contains all possible states of the system at the time of the 
initiation of control, is a closed subset of 6. The state equation consists 
of the kinematic equation and the dynamic equation. There are two types of 
perturbation: y^ represents target angular velocity; y^ represents distur- 

bances entering through the dynamic equation. The system will be said to be 
under the action of the nominal control law when y^ = 0. In general, the 
system may operate in any one of several modes. A system with several star 
trackers is typical. In certain mechanizations the form of control law 
depends on which star trackers are active. Hence, in such systems there are 
as many modes as star tracker combinations. The mode is represented in table 1 
by the mode variable y. It will be assumed that the number of modes is 
finite, that with each mode there is associated an open subset 0„ of 9, and 
that these subsets cover 0 . The perturbations y^ and y ^ are given in terms 
of intensities nj and n2 and forcing functions ui and U2 . The combined 
forcing function 11 = (iii ,02) is restricted to the set of piecewise continuous 
vector functions of time with values in U. This set is denoted in table 1 
by IJ. Finally, the functions f, n^, and n2 are assumed to be continuous on 
the closure of 0^ for each y. 

The following examples illustrate how specific cases may be described by 
the model of the type shown in table 1. 

Example 1- Suppose that initial attitudes of both the spacecraft and 
target are arbitrary, that initial spacecraft angular velocity is spherically 
bounded by a given constant (^amax> but otherwise arbitrary, and that target 
angular velocity is spherically bounded by m^max for all t > 0, but other- 
wise arbitrary. Then, the set of admissible initial conditions 0 q in 
table 1 may be defined by (see eq. (16)) 


0Q = {x: 0 < <|) ^ m, II 0)^11 < co^max^ 

and the target angular velocity perturbation may be accounted for by 


“^dmax* llujil <_ 1 


The form of the dynamic equation depends on the way torque is generated. 
Suppose that the spacecraft is being controlled by means of reaction wheels 
(see eq. (24)), and that the total angular momentum of the system hg = 0. 
Suppose, also that the feedback linking the sensors to wheel motors is given 
by 
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= z(R,io^,u) 

where z is a specific function. Then the function f in the dynamic 
equation in table 1 is given by 

f = Ja^z(R,Wa,y) + 

while the intensity of perturbation in this case is 

n2 = 0 

Example 2- Suppose that the situation is as in example 1, except that 
the control is by means of control moment gyros (see eq. (25)). Let the input 
e to the processor be given by 


e = z(R,Ua,u) 

and let the perturbation function be defined as 

n2 = max III - J"^hq(q)F(q) J^ll il z(R,o3^,p)ll , II U 2 II 1 (26) 

This gives the maximum deviation from the nominal control law z(R,a)^,lj). 

Then the system can be modeled by setting 

f = z + y^ 

The resulting model in table 1 represents the system in the sense that all 
possible motions of the systems are included in the set of all possible 
motions of the model. 

Example 3- Suppose that the spacecraft is controlled by means of an 
angular momentum exchange and storage device as in example 1, and that an 
angular momentum dumping scheme is employed. The corresponding dynamic equa- 
tion is given by equation (22) . Let the dumping torque be spherically bounded 
by and let the total angular momentum hg be spherically bounded by 

^smax- Then, assuming that (h^)‘ =-z (R,cj^, p) , the system can be represented in 
the form of table 1 by setting 

f = J"^z(R,o)a,M) + 


y I = n2U2 


13 



where ri2 is the 3x6 matrix 


Ti2 


[^max'^a 


"smax’ 


Ja S(o)a)] 


(27) 


and where the forcing vector is six-dimensional. 


llu^ll ^ 1 for i = 1, 2 


DESIGN OF CONTROL LAWS 



The direct synthesis of control laws by means of optimal control theory 
is often impractical because of the attending computational difficulties and 
because the resulting laws are difficult to implement. An alternative 
approach is to pick a feedback structure which on physical grounds is likely 
to result in adequate system performance, and then perform a global analysis 
to determine whether the performance is, in fact, adequate. The latter 
approach is taken in the present report, A method for generating a class of 
candidate control laws is presented in the present section. The system 
governed by any member of this class is shown to be asymptotically stable. 
Further global properties of the system may be determined by the techniques 
developed in the succeeding sections. 

Suppose that the target attitude is constant = £) so that the 
kinematic equation (13) is 


R = S(o)^)R 


(28a) 


and that the controls are adjusted so that the dynamic equation (22) is 


(La = z(x) 


(28b) 


where z is a control law to be selected, and x = (R,(ji)g^) is the state. The 
method for generating candidate control laws to be now discussed is based on 
the following simple example. 

Choose the magnitude of attitude error to be the following function of R 

m(R) = I (29a) 


14 



where the error angle 4> is defined by equation (16) . The time rate of 
change of m along any trajectory of the system (28) is, according to 
equation (18a), 


m = (4^c)^a)g^ (29b) 

Let the control law be given by 

z(x) = -())C - (29c) 

Then the system (28) is asymptotically stable on the set 


0 = 
0 


x: 


m 


1 t 
+ — ^ ^ 

2 a a 


< a 


for any a in the half-open interval [0, (1/2 ) tt^). The function 

V = m + (1/2 )(jo^cj is a Lyapunov function for the process: V is a positive 
a a * t t t 

definite, and its time derivative along any trajectory, V = <|)C 

is negative for ^0, and if cog = 0, then z ^ 0, unless also ip, or 
equivalently m, is zero. 

The control law (29c) has been generated from the error function (29a) in 
the following sense. The procedure was to pick a function m that character- 
izes the magnitude of the three-axis attitude error. Its time derivative 
turned out to be a linear function of Wg, namely ((J)c)^a)g. This function 
(covector) ((fic)^ was converted into a vector <tic and taken as the attitude 
error feedback portion of the control law. Finally, rate feedback was added 
for damping. 

Such construction can be carried out in general because the kinematic 
equation (28a) is linear in Wg. Thus, let s(R) be a differentiable function 
of R on the region of interest. If s is bounded for all s of interest, 
then becaiise of linearity of (28a), there is a matrix F(s) such that along 
any trajectory of the system s = F(s)o)a,- Choose a differentiable nesting 
function m(s) that assigns the notion of magnitude to attitude error in terms 
of s; that is, if a^b then {x: m[s(R)] ^ a) c {x: m[s(R)] .Sb). The 

time rate of change of m along any trajectory is m = m 5 (s)F(s)ma. If for 
all m ^ mo, ms(s)F(s) = 0 only at m = 0, then the system (28) controlled by 


z(x) = -F^(s)m^(s) - G(s,o)g)a)g 


(30) 


where G is a positive definite matrix, is asymptotically stable with respect 
to m on the set 0o = (x: m + (l/2)m|ojg j< mo>. A Liapunov function is 
V = m + (l/2)w^a)g. The time derivative, V = -a)|Go)a, is by assumption negative 
for (jjg ^ 0. For ug = 0, z / 0 unless attitude error feedback is also zero. 


15 



which by assiimption occurs only at m = 0. Hence, any initial condition in 
®o will be driven to (m = 0, = 0) . This may be, depending on m, a set in 

6 - see example 5. 

-(j) 

Example 4 - eigenvector control- Let m(R) = I g(q)dq, where qg(q) > 0 

•^o 

for q 0. This error function generates the following control law. 


z(x) = -g(<|))c - G(x)wa 


(31) 


The attitude error feedback vector is along the eigenvector c of the error 
matrix R, In particular, let g(<t>) = 2k sin[(l/2) (j>] , where k is a constant. 
Then the attitude error feedback becomes -ke, where e is the Euler vector 
defined by (19b) . 

Example 5- The problem is to_aline a unit vector b fixed in the 
spacecraft with the unit vector d fixed in inertial space. Let the space- 
craft coordinates of these vectors be denoted by b^ and d^, and let their 
inertial coordinates be denoted by bg and dg . In this case, ba and ds are 
constant and da is measured. A solution to this problem may be obtained by 
choosing the error function to be, for some positive constant kj, 

m = ki(l - d^b^) (32) 

Then, since da = Aasdg, and using the kinematic equation (10), it follows 
that 


lii = kid|s(o)a)b^ = -kid^S(b^)iOa 

If G in equation (30) is taken to be another positive constant k 2 , the 
following control law results 

z(x) = -kiS(ba)da - k 20 J^ (33) 

It may ^e noted that kiS(ba) is a constant matrix. The system will aline F 
along d for any initial condition in the set defined by m+ (l/2)u)ga)^ < k^ . 
Note that m is semi definite since rotation about d is irrelevant; hence, 
m = 0 on a set in 9. The chosen magnitude fun£tion can be written as 
m = 1 - cos 0, where 9 is the angle between b and d. So, m measures the 
angle between the t£o vectors . Suppose that m were defined as the magnitude 
of the difference b - d; that is, let raj = (1/2) (d^ - b^)^(da - b^) . 
Expanding this product and noting that b and d are unit vectors, it follows 
that mj = 1 - d^b^ = m, and the same cross product control results. 
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In certain cases linear control is undesirable because it may require 
excessive torques for large errors. Nonlinear gains may be introduced to 
account for torque saturation as follows. Let s = da - b^, and choose the 
magnitude function of s to be 


m 


= ^ ^ g^Cqidq 

1 0 1 


with qg^(q) > 0 for q 0. Then the time rate of change of m is 


m 


= S gi(Si)s. 


But, s = -S(da)tOa, which leads to the following control law. 


'g^Cdai - bai) 

z = S(da)| g2(da2 " ba2) ) - 

vS3(da3 “ ba3), 


In particular, if all gj^ are the same saturation function ksat(q,qg) which 
is linear for |q| qg, and saturates at the value k for |q| ^.qg, then the 
attitude error portion of the feedback is bounded on each axis by k. 


Example 6- It is assumed that complete attitude control is desired. Two 

unit vectors b^,b^ fixed in the spacecraft are to be alined with two unit 

vectors d^,d^. inertial space. The spacecraft is on target (i.e., 

R = I) when b^ = d^, i = 1, 2. It is assumed that the body coordinates of 

di are measured by on-board sensors (i.e., a sun sensor and a magnetic field 

sensor). Let s^ = d^ - b^, i = 1, 2, and choose the magnitude function of 

a a 

attitude error to be 


m = 


E 


1 



g^(q)dq + 


2 

J hi(q)dq 


with qg. (q) > 0 and qh^(q) > 0 for q 0. This error function generates the 
following control law 


z 


( giCdal 
S3^a3 


) 1 + 
a2^ I 

I=a3)/ 


/hlC4l 

S(dJ)| hjCdL 

VhsCdL 



Goja 
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Figure 3.- Arrangement of star trackers. Body axes are 
numbered; d* are lines of sight; are inner gimbal 
angles; 7^ are outer gimbal angles. 


If the attitude error feedback is 
zero only at the point R = I of the 
region m v, then the system is 
asymptotically stable on the region 

m + (l/2)o)ta)^ < V. 

Example 7- It is assumed that 
spacecraft attitude is measured with 
a set of star trackers whose arrange- 
ment in the body is shown in 
figure 3. The body coordinates of 
lines of sight are. 



d: = 



d. = 



The corresponding kinematic equations are. 



SYi 


-cYite^ 


0 

-1 




SY2 

-CY 2 t 62 





a 



/ ^^3 
\sY3tB3 


0 

1 




Consider, first, the use of all three trackers. Let s 
angle error, si = - g?, i = 1, 2, 3. Superscript o 
values. Choose 


be the inner gimbal 
denotes the target 



(q) dq 


where qgj^(q) > 0 for q ^ 0. This error function generates the following 
control law. 
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z 


(34) 



CY2 

SY2 

0 



/giC^i 

( 82(^2 - 
VggC^s 


^;a 

B°) |- 

6°)/ 


Gu)a 


If in the region m(s) ^ v, Y2 ■*■ "*'2 ^ ±^/2, then the system is 

asymptotically stable with respect to R = I, in the region 

m(s) + (l/2)u)t(i) < V. It may be noted that the resulting control law may be 

simple to implement: the elements o£ the gain matrix are provided by 
resolvers attached to outer gimbals of trackers. 

Now, consider only trackers 1 and 3. Let 

st = (6^ - bJ, Yi - yJ, Bg - 33, Y3 - Y°) 


and again choose the error function to be 


m(s) 


4 


E 


r^i 

J g^Cq^dq 


It generates the control law. 



(35) 


If in the region m(s) <. v, ±tt/2, 33 ±if/2, and Y^ + Y3 ±'’f/2, then 

the resulting system is asymptotically stable with respect to m on the set 
of initial conditions m(s) + (l/2)o)3u)^ ;! v. If, in addition the two lines of 

sight are independent, then m 0 implies R I . This control law is 
harder to implement than the preceding one because of the presence of the 
tangents of inner gimbal angles in the gain matrix. If the terms involving 
tangents are set to zero, the attitude error feedback becomes that used in an 
actual Orbiting Astronomical Observatory (OAO) . The effect on system perfor- 
mance of this change may be investigated by the techniques presented in the 
following sections (see example 10). 


The final example of this section illustrates the design using actual 
hardware. 
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Example 8- It is assumed that the spacecraft is similar to an OAO. The 
control is by means of reaction wheels. The angular momentum storage capacity 

of each wheel is h^^^x = 4.68 (N-m-sec) ; the torque capacity of each wheel 

motor is Tmax = 0.231 (N-m) . The moment of inertia of the main body is diag- 
onal: Jg = diag(999;1110; 1410) (kg-m^) . An angular momentum dumping scheme 

is assumed which maintains the total angular momentum spherically bounded by 
^smax “ 1.00 (N-m-sec). Hence, the wheels will not saturate if body angular 
velocity cog is kept spherically bounded by 


‘^amax ~ Ch„,gx - hg|^ax)/jj^^j^ 
= 2.61 (mrad/sec) 


where is the maximum principal moment of inertia of the main body. 

Finally, it is assumed that on-board sensors measure the attitude error R 
and body angular velocity Wg. The problem is to design a control law for 
stable attitude regulation. 

Choose the magnitude function of R to be 


ni(R) = J kisat(q,q )dq 
0 ^ 


It generates the control law 


z(x) = -kisat((|),(t)g)c - k20)^ 


(36) 


Let kj = Then llwgll'^O for II 2. (Ugj^gx. Hence, the control will 

not drive the wheels into saturation. Let k^ + = Tmgx/jjj^^^. Then if 

the control law is implemented by setting the wheel motor torques 
T™(x) = J^z(x), the motors will not be driven into saturation. Finally, 
choose the saturation point (jig so as to have damping c = 0.5 near 4> = 0 . 
This specifies the control law completely. It may be noted that since 


sat((f),(t)g)c 


sat(.^,(t)g) 
2 sin i() 




the control law (36) is well-behaved at (|) = 0. 
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Components of The control torque, 10 J Nm 



Figure 4.- A particular response of the system in 
example 8 in terms of magnitudes of torque, 
angular velocity and attitude error. 



Figure 5. - The response of the system in example 8 in 
terms of the components of the control torque. 


The response of the system to 
the initial condition, wheels locked 
prior to t = 0; $(0) = 2 (rad); 
ct(0) = 

“ 0.5(1,-1,-1) (mrad/sec) is 
plotted in figures 4-6. These curves 
were obtained on a digital computer. 

The response may be divided 
roughly into three parts. For 
0 ^ t ^ 100 seconds, the control 
generates pulselike torque to bring 
the spacecraft to its maximum angular 
velocity. For 100 ^ t ^ 750 seconds, 
the vehicle coasts with maximum veloc- 
ity toward the target; the small 
torques in this time interval proba- 
bly counteract the gyroscopic term in 
equation (24). For t >. 750 seconds, 
the control again generates pulselike 
torque to stop the spacecraft on tar- 
get. At the end of the transient, 
the angular momentum of the system, 
which was initially in the main body. 



t, sec 


Figure 6.- The response of the system in example 8 in 
terms of the components of the relative momentum 
of the reaction wheels. 
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This system was, further, breadboarded on an airbearing platform using 
reaction wheels, motors, and star trackers (see ref. 8). This simulation 
showed the design to be practical. However, before the design can be con- 
sidered to be complete, the following type of questions must be resolved. 

(1) How does the system respond to any admissible initial condition, and how 
well does it follow a moving target? (2) How significant is gyroscopic 
coupling? How sensitive is the system to (3) external torque disturbances, 

(4) variations in system parameters, and (5) changes in the form of the con- 
trol law caused by partial failures in sensors and torquers? Such questions 
may be resolved by means of the techniques presented in the following sections. 


GLOBAL ANALYSIS OF SYSTEM PERFORMANCE 


Suppose now that a single-mode attitude control system has been designed. 
(See Applications, Case 2, for the discussion of a multimode design.) The 
problem is to determine whether the proposed design will perform adequately in 
all possible control situations. A control situation is defined by an initial 
condition, Xq, and a pair of forcing functions, u = (ui,U 2 ); Xg is the state 
of the system at the time the control is initiated; u^ generates a motion of 
target attitude; and U 2 generates a time history o"f disturbance. The latter 
may be due to external torque, or, as in a sensitivity analysis, it may repre- 
sent variations in system parameters. The environment in which the system 
will operate is characterized by the set of admissible initial conditions 6 
perturbation functions ni and n 2 (see table 1) and the set of admissible 
forcing functions U. Every motion of the system, say x, is the solution of 
the state equation x = g[t,x,u(t)] for some Xq in 0q, and some combined 
forcing function ^ in IJ. That is, the state equation induces a transforma- 
tion of the Cartesian product set 0gxU onto the set X of all possible 
motions of the system. The elements of X will be denoted by x=x(t,x ,u) 
foi" t ^ 0. Global analysis is concerned with the properties of this set 
X. The problem is to characterize X, and then to devise an effective 
procedure for computing the chosen characteristic. 


Response Envelope 

The primary purpose of an attitude control system is to maintain the 
spacecraft on target attitude. Hence, that property of X is of interest 
which characterizes the overall behavior of attitude error. To decide at any 
instant of time how near the actual attitude is to the desired attitude, it is 
necessary to have a notion of magnitude of attitude error. This may be intro- 
duced by means of the attitude error function ra(R) as discussed in the pre- 
vious section. Although this function may be the same one that was used to 
design the control law, in general, it may be different. Each motion x of 
the system generates a corresponding time history m of the magnitude of 
attitude error. It is a positive scalar function of time which will be 
denoted by m(t,XQ,^) for t >. 0. Let M be the set of all such curves gen- 
erated by X, or equivalently, by the set of all possible control situations 
0gXU. The system response envelope, denoted by m**, will be defined as the 
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function of time which at each instant is the maximum of all values of m in 
M at that time; that is, for all t >: 0, 

m**(t) = max max m(t,Xo,u) (37) 

Xo^Qo 

Thus, the response envelope is a global property of the system. It indicates 
system responsiveness. For any admissible initial condition, target motion, 
and disturbance the attitude error at any time t >. 0 will not be greater 
than m**(t). 


State Space Interpretation of the Response Envelope 

The response envelope can be given the following state space 
interpretation. Let 0 and 0 q be represented schematically as in figure 7. 
Since m is independent of Wa. the surfaces of constant m(R) are nested 
cylinders in 0. The motion starting at a.nd forced by is a 

trajectory in 0. The cylinder being crossed at time t determines the 
associated value of m at t. 

The same initial state Xq but a different forcing function u will 
generate a different trajectory. Consider the bundle of all such trajectories 
emanating from Xq and generated by IJ. This bundle defines at each t >. 0 
a moving set of states that are reachable from Xq at time t. Such a set is 
shown crosshatched in figure 7. The maximum cylinder intersecting this set at 
t gives the value of m corresponding to the inner maximization in 
equation (37). In figure 7 it is denoted by m 2 . 

Now, consider the union of all such moving sets generated by all Xq^Qq- 
This union defines a moving cloud of points 0(t,0Q) shown schematically for 
t = 0 and t = ti in figure 8. Initially, the cloud coincides with Og. 


Boundary of the 
set of admissible 
initial conditions 



(Surface of constant m) 


Figure 7.- Motion of set reachable from Xq € 6 q . 
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Thereafter it moves in response to system dynamics. At each t >. 0 the 
maximum m-cylinder intersecting the cloud gives the value of the response 
envelope m^** at that time. In other words, the knowledge of the motion of 
the boundary of the cloud is sufficient for the computation of the respons*^ 
envelope . 


As an illustration of the preceding discussion consider the 
unforced linear second-order system 



Figure 9.— Motion of the cloud of states for the example. 



following 



with 0Q being a square centered at 
the origin and bounded by xi = ±1, 
X 2 = ±1, and choose arbitrarily the 
error function ra(x) = lxi|. The 
corresponding motion of the cloud is 
sketched in figure 9. The corre- 
sponding response envelope is given 
in figure 10. In this simple case 
the response envelope is determined 
by the motion of vertices A and B. 

Returning to the general case , 
consider the evolution of a piece of 
the boiondary of the cloud, which at 
time ti is shown schematically as 
ABC in figure 11. At a later time 
t 2 it transforms into A'B'B"C, 
which is the envelope of all moving 


B’ 



Figure 10.- Response envelope of system in example. Figure 1 1 . The Huygen’s construction of the boundary 

of cloud. 
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sets emanating from ABC. The image of ABC at some still later time t3 is 
given by the envelope of all moving sets emanating from A'B'B”C', and so on 
for later times. It can be seen that this is the standard construction of 
wave fronts. Thus, for every initial condition Xq on the boundary of Qq, 
there is a u with values on the boundary of U(t) for every t ^0, such 
that the resulting trajectory remains on the boundary of the moving cloud 
for t ^ 0. Note that the system was assumed to be in a fixed mode. 

Suppose there is a function V(t,x) such that VCt,x) >0 outside, 

V(t,xj = 0 on the boundary, and V(t,x) < 0 inside the cloud. Consider the 
differentiable portions of the boundary and let V(t,x,uj denote the time rate 
of change of V along a trajectory. Note that since V depends on trajec- 
tory, it depends on u. According to the preceding discussion, the boundary 
is characterized by two properties: (1) it is part of the cloud, that is, at 

each point on the boundary there is a such that V(t,x,u) = 0, and (2) no 

trajectory can penetrate the boundary outward, that is, at each point on the 
boundaiy and every uGU, V(t,x,uj ^ 0. Therefore, V satisfies the following 
equation on the boundary 


max V(t,x,u) = 0 (38) 

uGU 

But, V = + VxX, and x = g(t,x,u) where the subscripts indicate partial 

differentiation. Hence, equation (38) can be expressed as follows: 

Vt + H(t,x,Vx) = 0 (39) 

where 

H(t,x,V^) = max Vx(t,x) g(t,x,u) (40) 

uSU 

Thus, the differentiable portions of the boundary satisfy the Hamilton-Jacobi 
equation (39) with the boundary condition (x:V(0,x) ^ 0} = 0^. The correspond- 
ing response envelope is given for t >. 0 by 


m**(t) = max m(x) (41) 

(x:V(t,x) = 0} 

That is, m**(t) is the maximum cylinder intersecting the cloud at t. 


COMPUTATION OF RESPONSE ENVELOPES 


According to the preceding discussion, the computation of a response 
envelope involves the solution of the Hamilton-Jacobi equation (39) . The 
response envelope is given in terms of this solution by equation (41) . Two 
methods for solving (39) are discussed in the present section. One, an 


25 



approximate method, is based on the Liapunov theory of stability. The other, 
an exact method, is based on the theory of optimal control. 


Approximate Computation of Response Envelopes 

Suppose V^(t,x) = 0 is a smooth surface which at each t ^ 0 encloses 
the moving cloud without necessarily being its boundary; that is, suppose that 
for each t >. 0, 

{x:V(t,x) <. 0} c {x:V*'(t,x) ^ 0} 

Then property (1) of the exact boundary may be dropped with the result that 
V'*’ satisfies the following inequality which is characteristic of Liapunov 
functions 


+ H(t,x,Vp i 0 (42) 

This is the Hamilton-Jacobi equation with equality replaced by (^ . The 
boundary condition is 

(x:V'^(0,x) .s. 0} 0Q (43) 

Equation (42) must hold for all t >. 0, and all x in 6 such that 
V^(t,x) = 0. 

It is very easy to construct solutions to (42) . Simply, let 

V^(t,x) = Vi(t,x) - V 2 (t) (44) 

where Vj is such that for some finite a. 


{x:Vi(0,x) i a) 3 


(45) 


and where V 2 is the solution of the following ordinary, first order, scalar 
differential equation with V 2 ( 0 ) = a. 

V 2 = max W (46) 

{xee:Vi = V2} 

where 


W = + H(t,x,Vj^) (47) 
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Then so defined solves (42) as can be seen by direct substitution. The 

corresponding approximate response envelope is given by 

m'^(t) = max ra(x) (48) 

{x60:Vi = V 2 } 

By construction, m^(t) >. m**(t) for all t >. 0 . For this reason such an 
approximation will be called an upper estimate of the response envelope. Such 
an estimate is useful because it may serve as a basis for accepting a proposed 
design: under no circumstances can the attitude error be larger than m'^(t) 

at any t ^ 0. Any function satisfying the boundary condition (45) may be 
vised to confute an upper estimate. Of course, the fidelity with which m 
represents m** depends on the choice of Vi. A poor choice will result in 
an overly pessimistic estimate of system performance. The following simple 
example illustrates the above discussion. The selection of Vj for attitude 
control systems is discussed later (see eq. (51)). 

Consider the second-order system 


xi = -xi 


X 2 = (2e - 2 )x 2 + X 2 U 

where the forcing function |u| ^ 1. Suppose that the set of admissible 
initial conditions is the unit squap, namely, 6p= (x: Ixi] <. 1, IX 2 I ^ 1), 
and that the error function m = (Xj^ + X 2 ) 


2 ^ 1/2 


Let 


Vi = XI 


2 

X2 


Then condition (45) is satisfied for a = 2. The time derivative of Vi 
along any trajectory is 


Vi = -2x^ + 2(2e - 2 + u)x 2 


Hence , 

W = -2x5 2(2e‘^^ - l)x| 


and the maximum of W on the boundary Vi = V 2 is (4e 2)V2. Therefore, 

the differential equation (46) is 


V 2 = (4e - 2)V2 


whose solution for 


V 2 ( 0 ) = 2 is 


V 2 (t) = 2e 


2(l-t-e”^^) 
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So that the upper estimate is 


m 


■Ct) = 


On the other hand, in this simple example the exact response envelope can be 
obtained analytically as 


m 


I* * 


(t) = 


-2t 

e + e 



J 

Both curves are shown in figure 12 
for comparison. The difference 
between m** and m'*’ is due to the 
difference between the moving cloud 
and the approximating set Vi - V 2 ^ 0 . 
Thus, the cloud becomes squashed 
along xi; while, the approximating 
set remains circular. 

The calculations in the above 
example were sufficiently simple to 
be carried out by hand. For practi- 
cal systems, these calculations will 
most likely have to be done on a 
computer. An outline of a possible 
computer program is given in table 2. 


Figure 12. Response envelope and an upper estimate. 

of m"^ on the interval 


It is assumed that Vi 


and '^a” 


estimate 
estimate 
grid G 


satisfying (45) have been selected. 
The computation results in an 
_ 0 ^ t <. T, for some chosen T. The 

_ IS assumed to be sufficiently accurate when a refinement of the 
causes no significant changes in 


m 

in 


m. 


TABLE 2.- COMPUTATION OF UPPER ESTIMATES 

Step 1. Set V 2 CO) = a 

Step 2. Cover the surface Vi(tk,x) = V 2 (tk) by a grid G. 

Step 3. Compute W from equation (47) and m from its 
defining equation, and maximize both over G, 
denoting the maximum of m by m. 

Step 4. Step forward: VaCt^+i) = V 2 (t,^)W(tj^)At. 

Step 5. Repeat seeps 2 to 4 for all tj^ in [0,T]. The 
result is an estimate m of m'*' over [0,T], 


It may be noted that table 2 represents the computational procedure that 
IS followed in a typical stability analysis by means of Liapunov functions. 
The only modification, aside from the maximization of m, is that here the 
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Liapunov conditions are tested at Vi = V 2 with V 2 a computed function of 
time, rather than Vi = for some preassigned collection of numbers {v^}. 


Exact Computation of the Response Envelope 

Consider now the exact solution of the Hamilton- Jacobi equation (39) . 
Suppose that V(t,x) is differentiable for all x of interest and all t >. 0. 
Let the trajectory x(t) = x(t,XQ,u) lie on the boundary V(t,x) = 0, that is, 
V[t,x(t)] = 0 for all t ^ 0. Let p(t) be the outer normal to the moving 
surface along x(t) , that is, p(t) = V^[t,x(t)] for all t ^ 0. The differ- 
ential equation satisfied by p may be obtained as follows. Let 
x*(t) = x(t,x*,u*) be a neighboring trajectory, which also lies on the moving 
boundary, and°let 6x(t) = x*(t) - x(t) . For sufficiently small 5x 
V(t,x*) = V(t,x) + p^Sx. Hence, pt(t)6x(t) = 0, for all t i 0. Consequently, 

p^6x + p^(6x)' = 0 

But the system state equation is x = g(t,x,u). Hence, 

(6x) ■ = g^(t)6x + gu(t)6u 

where the coefficient matrices are evaluated along x(t) , and u(t) , and 6u 
is such that p^g^^u = 0 because the neighboring trajectory x* is also on 
the moving boundary. Therefore, 

[pt + ptg^(t)]6x = 0 

since 6x is an arbitrary vector in the tangent space of V(t,x) =0 at 
X, p + g^P = kp for any scalar k. Choosing k = 0, one obtains the follow- 
ing differential equation for the normal p 

P = -gx(t)P 

Thus, the motion of a planar element which is given at t = 0 by position 
Xq and normal p^ satisfying V(0,Xo) = 0 and p^ = V^CO ,Xq) , respectively, 
satisfies the following standard equations of optimal control, which are 
derived rigorously in reference 9: 


X = g(t,X,u) 

(49a) 

p = -gj(t,x,u)p 

(49b) 

u = argmax H(t,x,p) 
uG u 

(49c) 
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where H = p^g(t,x,u). Equation (49c) expresses the Huygen's construction 
(fig. 12) of the boundary. The time history of attitude error m* corre- 
sponding to this planar element is thus a function of only x^ on the initial 
boundary. For each t >. 0 the response envelope is given by 


m**(t) = max m*(t,Xo) (50) 

(xo:V(0,Xq) = 0} 

This coiq)utational procedure is outlined in table 3. It is assumed that 
Sq “ (x:V(0,x) ^ 0}, that V(0,x) is given, and that it is smooth. The 
estimate is assumed to be sufficiently accurate when a refinement of the grid 
G causes no significant change in m. Pontryagin's maximum principle 
(ref. 9) guarantees that m will converge to the exact response envelope m** 
with grid refinement if equations (49) have unique solutions for given initFal 
conditions. This will be the case if, for example, in table 1, g is differ- 
entiable, perturbation functions n^ and n£ are smooth and have maximal rank 
almost everywhere, f in the dyneimic equation is invertible with respect to 
y^, and U = {urilull <_ 1}, 


TABLE 3.- COMPUTATION OF THE RESPONSE ENVELOPE 


Step 

1. 

Cover the initial surface V(0,x) =0 by a 

grid G. 

Step 

2. 

Set m = 0 on the time interval [0,T], 


Step 

3. 

Set initial conditions (x,p) = [x^ ,V^(0,X4 ) ] 

i e G. X r 

for 

Step 4. 

Solve^the canonical equations (49) storing 
max[mj^(t) ,m(t) ] . 


Step 

5. 

Repeat steps 3 and 4 until G is exhausted, 
resulting m is an estimate of m** on the 
computation interval. ~ 

The 


It may be noted that this computation of the response envelope is 
essentially no more complicated than the maximization of m(T) for a single 
initial condition Xq, In the latter case all directions of initial p must 
be tested. This is equivalent to 6^ being an infinitesimal sphere about 
the initial state Xq. 

The main advantages of the computation outlined in table 3 are that 
there is no need to guess a Liapunov function, and that, at least in the 
limit, the exact response envelope is being computed. On the other hand, the 
con^jutation outlined in table 2 does not require repetitive solution of the 
system state equation or its adjoint. Hence, the conditions on the state 
equation are much weaker. Of course, there is also the advantage that in 
certain nontrivial cases the required computation can be carried out by hand. 
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If the computer is to be used, the computation time required is of 
practical interest. This can be estimated by assuming that the set of admis- 
sible initial conditions 6 q is a sphere. Assuming that there are sub- 

divisions of each coordinate interval, there are Ng = 2nN^“ grid points in 
G. This is the number of computations involved in step 3 of table 2. 

Assuming that there are N-t points along the time interval, the total number 
of computations is N = 2nN^”^N^. On the other hand, in table 3 there are Nq 
time histories, each requiring N-^ computations. Thus, in either case there 
are N = 2nN"“ computations. For example, if Nx = 10, = 100, then 

N = 1.2x10^ if n = 6, and N = 6xl0‘* if n = 3. Assuming 100 microseconds 
per computation, the computation time is of the order of 3 hours for n = 6, 
and only 10 seconds for n = 3. The large reduction in computation time 
accompanying the reduction in the dimension of the state space motivates the 
following discussion of comparison models. It will be seen that the perfor- 
mance of an attitude control system can be compared with that of a spherically 
symmetric model whose state space is essentially three-dimensional. 


Comparison Models for Attitude Control Systems 

From the practical point of view, it is very desirable to be able to 
trade accuracy for reduced computer time in a meaningful way. A useful 
approximation to the response envelope is an upper estimate m'^ such that for 
all t 0, 


m+(t) ^m**(t) 

As noted previously, such an estimate may be used as a basis for accepting a 
proposed design: for no combination of possible initial condition, target 

motion, and disturbance will the magnitude of attitude error be greater than 
xj^(jtQ 3 .ted by the upper estimate. Such an estimate may be computed using a 
comparison model having two properties. First, the comparison model must be 
sufficiently simple that the computation of its response envelope is practical. 
Second, it must be known analytically that this response envelope is an upper 
estimate on the response envelope of the given system. A way for constructing 
comparison models will now be discussed. 

Suppose that two initial states Xi and X 2 happen to be such that for 
every admissible forcing function there is an admissible forcing function 

u^ such that for all t ^ 0, 


m(t,xi,u^) = m(t,X 2 ,u^) 


and conversely. Then one may consider the two states to be equivalent for 
the computation of the response envelope. It would be a waste of computer 
time to include both states in the grid G. For efficient use of the computer, 
the grid should consist of only the representative states, each representing 
its equivalence class. For example, suppose that the state space is 
n-dimensional, that the state equation is 

X = -X + u 
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that the set of admissible initial conditions is 


6q = {x: x^x -2^0} 


that for all t >. 0 , 


U = {u: Hull ^ 1} 

and that the magnitude function is 


m(x) = II xll 


Then it is sufficient to include only one point in G, say x^ =(i/2,0,. . . ,0) 

This simple example suggests that the concept of state equivalence is 
potentially useful for speeding the computation of response envelopes. Its 
actual usefulness depends on the ease with which equivalence can be determined 
The equivalence of two states can always be determined during the computation 
of the response envelope. Of course, this is not very helpful. Efficiency 
is achieved only if equivalence is determined before the computation is ini- 
tiated. It seems that for an arbitrary system such an a priori determination 
of equivalence is difficult. But, consider the situation from the other end. 
That is, start by choosing a partition of the state space, and construct a 
model whose set of equivalence classes coincides with this partition. Then 
adjust this model so that the set ^ of its possible motions includes the 
set )( of all possible motions of the given system. Then the response 
envelope of the model can be computed efficiently, and it will be an upper 
estimate of the response envelope of the given system. The desired trade-off 
between computer time and accuracy is thus accomplished. A fine partition 
will result in small saving of time, but the estimate will be close to the 
response envelope of the given system. In fact, for the finest partition, 
namely identity, no time is saved, and no error is made. As the partition is 
made coarser, equivalence classes become larger, computation time smaller, and 
the estimate more conservative. Of course, if the con^arison model happens to 
be the exact model of the given system, time is saved without loss in accuracy, 

A convenient way to define a partition is by means of a group of 
transformations. In that case two states xj and X£ are equivalent if there 
is a transformation taking xi into X 2 . A partition is obtained because a 
group has an identity (reflexivity) , an inverse Csymmetry) , and closure 
(transitivity) . 

In summary, a comparison model may be constructed for a given system as 
follows. Based on physical insight, choose a group of transformations. The 
choice defines a partition on the state space. Construct a model whose 
equivalence classes give this partition. Adjust the model so that the set of 
its motions )C^ includes all possible motions of the given system. The 
result is a comparison model of the given system. 

Now consider a comparison model for attitude control systems. The 
system state x can be represented by (e,u)g^) where c is the Euler vector 
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defined by equation (19b) and is body angular velocity. Then, x can be 

considered either as one 6-dimensional vector, or as a pair of 3-dimensional 
vectors. Consider the set of transformations t with elements where for 

each rotation matrix A, 


= (Ae,Ao)a) 

T is a group. Two states xi = and X 2 = equivalent if 

the triangle formed by and is congruent to the triangle formed by 

£2 and Hence, the only properties of the initial state that matter in 

the computation of the response envelope are the length of r, the length of 
03 ^, and the angle between these two vectors. Comparison models generated from 
T will be called spherically symmetric. One such model is given in table 4, 


TABLE 4.- MODEL OF SPHERICALLY SYMMETRIC SYSTEMS 


State space 

Region of operation 

Admissible initial conditions 

State equation 

Kinematic equation 

Dynamic equation 

Perturbation 

Target velocity 


X = X = Cc,u)a) 

0 = { X : 11 c II < 2 } 

Go = {x: IU-1P-+ VilUall^ < V 2 < 4) 

1 = t s(wa + 

Wa = - [f(ll E lOaic + a^wa] + 
yj = asUiCt) 


Disturb ance 

Admissible forcing functions 


Y 2 = f(ll ell) [-a 2 EUoCt) + a3S(E)u2(t) ] 

II Ui (t)ll 2 + juoCt) 1^ + II U2(t)ll 2^ 1, for all t2iO 


Magnitude 


of attitude 


error 


mfx) = (1) 


The state space is six-dimensional. The state is represented by the 
Euler vector e and body angular velocity ua- Points with Hell = 2 are 
excluded from the region of operation 6 because the kinematic equation (20b) 
is singular there. The set of admissible initial conditions Gq is an 
ellipsoid whose shape and size are determined by constant scalars Vi and V 2 . 
The angular acceleration is a sum of the nominal control law and a perturba- 
tion. The nominal control law is a weighted sum of the -Euler vector and 
angular velocity. The perturbation consists of two terms: one acts parallel 
to c; the other acts perpendicular to c. The forcing vector u is 
seven-dimensional. Three components uj are used to generate target velocity; 
one component Uq is used to generate perturbation parallel to e; and three 
components U 2 are used to generate perturbations perpendicular to e. The 
combined forcing vector is spherically bounded by 1. The intensity of per- 
turbations is determined by the constant scalars ao,a 3 , and as which are 
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assumed to be greater than zero. The magnitude of attitude error is the error 
angle (p. When the spacecraft is on target, <(» = 0; otherwise it is between 
0 and TT . 

That this model is spherically symmetric can be seen by considering the 
effect of any from i. Thus, t^( 6) = 0, and = 0^. The kinematic 

equation is spherically symmetric because 

(Ae)‘ = 

A£ = j AS(o)^ + y^)A^Ae + i nACoj^ + y^) 

= -|s[(Awa) + (Ayj)](Ae) +yn[(Awa) + (Ay^)] 

and llA/jll = llyjl = a 5 lluill. Spherical symmetry of the dynamic equation can be 
shown similarly. Finally, (p(ARA^) = (p(R) because the trace is an invariant 
under a rotation transformation. So, for any two initial states on the bound- 
ary of Gq, say xi and if there is a rotation matrix A such that 
^2 ~ > ^hen for any u^ there is a u^ such that for all t >: 0, 


4>(t,Xi,U^) = <J)(t,X2,u2) 


Consequently, the model is spherically symmetric, and the grid G in either 
table 2 or table 3 may consist of only representative states of the form 
X = where 

( cos 
sin 
0 

/Vo \1 /2 


0 ij; ^ 7T 

Thus, the computation of the system response envelope requires a maximization 
over only two parameters, w and ip. This is practical. 

The parameters aj, in table 4 were assumed constant in order to simplify 
discussion. It is clear that these parameters may be allowed to be functions 
of Hell, lltOg^ll, and as well as time. In addition, perturbations along 0)^ 

and perpendicular to it may be included. Thus, the condition of spherical 
symmetry is not so restrictive as might appear from table 4. If the given 
system is spherically symmetric, then the recognition of this fact can greatly 
speed the computation of its response envelope. If the given system is not 
spherically symmetric, then it can be represented by a spherically symmetric 



e - (V2 - v^w^) 


2 ^ 1/2 
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model by absorbing the asymmetry into perturbations. Then the response 
envelope of the comparison model will be an upper estimate of the response 
envelope of the given system (see example 10) . 

In the computation of upper estimates by means of Liapunov functions as 
outlined in table 2, it is necessary to give a Vi function. One such func- 
tion found to be useful in practice is given by the following equation: 

vi(x) = gi(<^)d^ + 2a^(l - cos ^ <1.) + ^ lluiall^ + J ^^‘"a 

where g((j.) = f[2 sin(l/2) 4>]/ [2 sin(l/2) ((>] , and f is the function appearing 
in the dynamic equation in table 4. It may be noted that this V| is spheri- 
cally symmetric: for any rotation matrix A, Vi[t/v(x)] = Vi(x). Its form may 

be thought of as a natural extension to three axes of the Liapunov function 
commonly used in the analysis of single axis servos. Thus, the first two 
terms on the right of equation (51) depend only on the magnitude of attitude 
error. The next term depends only on the magnitude of angular velocity. The 
last term represents coupling, which depends on these magnitudes and on the 
angle between the error axis and the angular velocity vector. 

Applications 

The following two examples illustrate the use of spherically symmetric 
models for the computation of response envelopes by means of procedures 
outlined in tables 2 and 3. 

Example 9- This example illustrates the computation of upper estimates of 
response envelopes by means of Liapunov functions (table 2) . Consider the 
system discussed in example 8, page 20. The nominal control law is given by 
equation (36). It is spherically symmetric. The time history of the error 
angle corresponding to a particular control situation is shown in figure 4. 

Now global behavior of this system will be considered. 

To simplify the discussion, let the system be normalized as follows: 
time t -► t/(jJamax; angular velocity 0 )^ • Wamax- dynamic 

equation in the absence of perturbation is 


<^a = - 7^ sat((t>,<|)s)c - ^ 


(52) 


The following upper estimates were computed using the Liapunov function, 


Vi(x) 


•'n 


sat(<t>,4s)d(t) + 


^1 - cos j <t>5llw3^1l^ + sin c'^^Ua (53) 


which is a special case of equation (51) . 
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Case I - single-mode nominal system- For this case the set of admissible 
initial conditions is assumed to be given by 


= {x: (}) i 2, II w II i 1} 


(54) 


It is also assumed that the dynamic equation is unperturbed (i.e., total 
angular momentum hg = 0 for all t >. 0) and that the angular velocity of the 
target is spherically bounded by a fraction b of the maximum angular veloc- 
ity allowed for the spacecraft (i.e., lU^H ^ for t >.0). 


are given in figure 13 
II (jOall > 1, points with 


Note from (52) and (54) that since llaial 


The results 
< 0 for 


1 


may be excluded from step 2, table 2. That 
is, one needs to consider only that 
part of the Liapunov surface which is 
inside the cylinder II a3a.ll — 1* "The 
curve b = 0 indicates the respon- 
siveness of the nominal system to 
step changes in target attitude. 

Thus, for any admissible initial con- 
dition, the attitude error will not 
be greater than the value indicated 
by this curve. It can be seen that 
the system is not only asymptotically 
stable on but it is essentially 

on target ((|) ^0.01 rad) no later 
than 3/a)a,max = 1150 seconds after the 
initiation of control. Curves with 
b > 0 indicate how well the system 
follows a time varying target. Thus, 
the curve b = 0.2, for example, 
shows that for any admissible initial 
condition and any target motion with 

angular velocity bounded by = 0.322 mrad/sec, the attitude error will 

not be greater than indicated by this curve. It is emphasized that a curve in 
figure 13 is not the response to a particular control situation. Rather, it 
is a global description of system behavior under all possible (there are 
infinitely many) control situations. 

Case 2 - multiple-mode nominal system- In this case the set of admissible 
initial conditions is assumed to be given by 



Figure 1 3.— Global response of nominal single-mode 
system for several bounds on target angular 
velocity. 


00 = (x: <j) £ IT, llwall 1} (55) 

It may be noted that this set includes points at which the error axis c is 
double valued, so that the control (52) is undefined there. However, the 
system may be controlled using three modes as follows. For a fixed 6 choose 
an "a" so that the maximum of 4) on the set 

01 = (x: Vi(x) :< a, llwall i D (56) 
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( 

is IT - 6 , and let the maximum of <j> on = 0 } be denoted by 41 ^^. 

If x(0) is in 01 , let the system be controlled by (52). Otherwise, apply 
maximum angular acceleration antiparallel to until either Oj is 

entered or ua ~ This mode lasts for at most II wg(0) Il/w 3 jj)a^x seconds. If 
01 is entered, let the system be controlled by (52). Otherwise, offset the 
error attitude R by A4> = rr - by introducing a fictitious change in 

target attitude. This brings x into 0i where the system can be controlled 
by (52) while the offset is removed with angular velocity which is bounded by, 
say. biu)..™.,.,. The effective target velocity in this mode is bounded by 
(b + bi)a)amax- The offset will be removed after (tt - <|>m) / Cb i^amax) seconds. 

Thereafter, the angular velocity is spherically bounded by ba)amax> as in 
case 1. The plots (fig. 14) show the behavior of the resulting system with 
d = 0.01, bi = 0.1, and w^max = 2/4>s = 20 rad/sec^. Note that curve b = 0 
shows the regulator is asymptotically stable for all attitude errors. 

Case Z - pevtuTbations •in the dynamia equation- In this case it is 
assumed that the dynamic equation is given by 

J)a = - sat(4),<t)s)c - ^^a ^2 

where the first two terms correspond to the nominal control law (52) and y^ 
is the perturbation. The set of admissible initial conditions is assumed to 
be given by ( 55 ), and the target attitude is assumed to be stationary. 


Figure 15 shows upper estimates due to perturbation of the form 



Figure 14.- Global response of nominal multimode system Figure 15.- Sensitivity to gyroscopic coupling, 
for several bounds on target angular velocity. 
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This perturbation is a (normalized) symmetric approximation of the 
gyroscopic term in equation (24) with llhsll £ bh^ax- The curve b = 0.3 shows 
that for the case considered (i.e., OAO) , gyroscopic coupling is not very sig- 
nificant even when the system is loaded with as much as 30 percent of its 
angular momentum storage capacity. 

Figure 16 shows the performance of the system with an angular momentum 
dumping scheme (see example 3, p. 13). It is assumed that the total external 
torque is spherically bounded by 0.1Tmax» and that the dumping scheme maintains 
the total angular momentum of the system spherically bounded by 0.3hniax- 

Figure 17 shows the sensitivity of the system to spherical errors in 
commanded acceleration. The perturbation is assumed to be given by 


y2 = b[llsat((t>,<|)g)c + a)all/<|)s]u2 


From this figure one may conclude that spherical errors of the order of 
10 percent affect the performance little. This means, for example, that 
10 percent changes in moment of inertia, motor and power amplifier gains, or a 
misalinement of the motor-wheel pairs with respect to the body axes of about 
3“ is not detrimental to system performance. Even when such errors are large 
enough to cause 30 percent error in acceleration, the system remains asymptot- 
ically stable. If the system were controlled by means of control moment 
gyros, the plots in figure 17 would indicate system sensitivity to partial 
failures in the gyro package. The corresponding b may be taken to be (see 
example 2, p. 13) 



Figure 1 6.- Sensitivity to gyroscopic coupling and Figure 1 7.- Sensitivity to spherical errors in acceleration 
external torque. 


38 




Figure 18.- Sensitivity to spherical errors in attitude 
error feedback. 


Figure 18 shows the effects of 
spherical errors in attitude error 
feedback. The perturbation is 
assumed to be given by 

= b[sat((|>,<t)s)/(t)g]u2 

The plots in this figure may be used 
to determine system sensitivity to 
errors and partial failures in the 
attitude sensor. 

Example 10- This example 
illustrates the use of spherically 
symmetric models and the procedure 
outlined in table 3 to compute upper 
estimates of response envelopes for 
systems which are not spherically 
symmetric . 


Consider the system discussed in example 7 (p. 18). Spacecraft attitude 
is measured with star trackers, and the difference between the actual and 
commanded gimbal angles is used for attitude error feedback. Let the control 
law be given by equation (35) in which the functions g^ represent hard 
saturation. In addition, let the terms in the gain matrix involving the tan- 
gents of inner gimbal angles be set to zero, and let the multiplication by 
rjixs matrix be followed by another hard saturation of each component. The 
resulting attitude error feedback is shown schematically in figure 19. Thus, 
the gimbal angle errors are clipped at 0.1 rad, passed through a pin matrix 
which is a function of the outer gimbal angles, and then again clippd. The 
result g(R) is the attitude error feedback. The dynamic equation is assumed 
to be the following linear combination of g(R) and body angular velocity 


= -lOg(R) 





Kinematics 


Star trackers 

Figure 19.- Attitude error feedback used in example. 


- 

The problem is to determine the 
behavior of the system on the set of 
admissible initial conditions given 
by 

0Q = {xrllell^ + llajall^ <. D 

The feedback g(R) is highly 
nonlinear, and it is not spherically 
symmetric. However, it can be repre 
sented by a spherically symmetric, 
smooth function with perturbations. 
Thus , in the range 0 <. <|) <. 1 rad , 
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g(R) - -10[1 + (10llel|)2] ^ ^ [e + O.SSeUq - 2.5S(e)u2] 


with Up + ||u2ll^ < 
computer.) Hence, 
32 = 5.5, 33 = 25, 
1.0 


1 . (This representation was determined on a digital 
in table 4 , vi = V2 = 1, f = 10[1 + (1011 e||)2]-i /2^ ^ ^ 

and 34 = 10. Figure 20 shows the corresponding response 

envelope computed by means of the 
procedure outlined in table 3. 

(Eq. ( 49 c) was made nonsingular 
almost everywhere by setting 
35 = 0,001 and rec 



Figure 20.- Global response of the system. 


II Uj 


+ ut 


II U2 II 


quiring that 
^ < 1 .) As 


can be 


seen from the figure, the system is 
asymptotically stable on 6 q, and it 
is essentially on target after three 
units of time for any admissible 
initial condition. 


CONCLUDING REMARKS 


approach to the design and global analysis of three-axis, large angle 
attitude control systems has been presented. The approach is general^in tL 

noUpnt special properties of particular^ystem com- 

ponents, but, rather, on properties common to all attitude control systems 
By making use of the well-known properties of three-dimensional ^otaUonr’ it 
was possible to apply the general techniques of control system theory to 
develop a practical desip and analysis procedure for such systems . Attitude 
rror, a kinematic equation, and a dynamic equation were formulated in a way 
that IS convenient for the study of attitude control systems, and were 
ec e in a pneral mathematical model of such systems. The notion of 
distance in attitude between spacecraft and target was introduced by means of 
attitude error functions. It was shown that such functions may be usTdTo 
generate asymptotically stable control laws. In addition, such functions may 
enveloped overall system behavior by means of response 


• interpretation of the response envelope was given, and the 

noted^^’T Liapunov's second method and optimal Lntrol theory was 

h.c r computing the response envelope were presented. 

One, based on Liapunov's method, is approximate and gives upper estimates on 
the response envelope. The primary advantage of this procedure is that few 
^ requirements are imposed on the system. The disadvantage is that 
there is no direct way to construct Liapunov functions. The second procedure 
based on the theory of optimal control, is exact and direct, but U imposer 
more conditions on system dynamics , ^ 

The computation time required by either procedure depends on the 
imension o t e state space. The concept of spherically symmetric comparison 
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models was introduced as a means for reducing the effective dimension of the 
state space from 6 to 3. This reduction results in a large saving of computer 
time. Any attitude control system with six-dimensional space can be compared 
with a spherically symmetric model by absorbing the asymmetry into perturba- 
tions. Of course, if the given system is strongly asymmetric, the upper 
estimate obtained will be overly conservative. 

The examples included in the report suggest that the proposed design and 
analysis technique is useful. 


Ames Research Center 

National Aeronautics and Space Administration 

Moffett Field, Calif., 94035, October 20, 1970 
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APPENDIX A 


METRIC PROPERTIES OF THE cp FUNCTION 


The (|>- function is defined for any rotation matrix R by 
(()(R] = arc cos I [trace (R) - I]1 

[0,^] f 

R may be interpreted as a rotation from d-basis into a-basis. 
paths from I to R. Each satisfies the differential equation 

R = S[oj(t)]R 

for some piecewise continuous w. In addition, R(0) = I and R(tr 
some fixed t£. It will now be shown that for all such u>, * 

rtf 

<))(R) < I lla)(t)lldt 

*^o 


The Hamiltonian is 


H = trace [pts(u))R] + pjla)!l 


and the adjoint equation is 


Po = 0 

P = S[o)(t)]P 


Thus, for any R and P have the same transition matrix '3>Ct). 
R(t) = 4>(t) , and P(t) = $(t)P^. Hence, for 0 t ^ tp 


trace [pts(w) R] = trace $] 

= trace [pts ] 

= 2o)t$k 


where k is a constant. Therefore, 

H = 2o)^<J>k + p II 0)11 
o 


Consider all 

) = R for 

(Al) 


That is, 
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and the optimum w is col inear with ^k, that is colinear with R(t)k. But 
this means that the direction of u) is fixed in the d-basis. Therefore, 
m(t) is at each t the eigenvector of R(t) , and the conclusion (Al) follows. 

The second property of ((> is the following. For any rotation matrices 
A and B, 


i <J)(A) + <t>(B) (A.2) 

Suppose the contrary, and denote AB^ by C and B by D^. Then it would be 
true that (KC) > ((>(0) + (|)(CD^). That is, the angle of the composite rotation; 
from I to D, followed D to C, is smaller than the angle of direct rotation 
from I to C. This, according to (Al) is impossible. Hence (A2) is true. 

Finally, consider the set of all rotation matrices. For any A and B in 
this set define 


<()(B,A) = c|)(AB^) 


(A3) 


The function 4>(B,A) so defined is a metric on the space of three-dimensional 
rotations. Indeed, (i) 4>(B,A) is positive; (ii) 4’(B,A) = 0 if and only if 
A = B; (iii) (|.(A,B) = 4>(B,A); (iv) <^(B,A) + <{-(A,C) i<j.(B,C). The triangle 
inequality holds because 

<()(B,C) = <()(CB^) = 4)[CA^(BA^)^] ^ <t>(CA^) + (KBA^) = (|>(A,C) + <|)(B,A) 
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APPENDIX B 


KINEMATIC EQUATION IN TERMS OF THE (4>,c) PARAMETERS 


From equation (16) it follows that 


-sin 4'4' = -j trace (R) 


But, according to (13), R = S(u))R. In addition, for any y in E^ and 


3x3 matrix A = (a. .) , 


trace [S(y)A] = -y^l ~ a 


a.23 - a32 


ai2 - a2i 


as can be checked by expanding both sides. Hence, 


( T23 - 1*32 
^31 - 1-13 
ri2 - T21 


which on using (17) gives 


= w^c 


To get (18b) note that Rc = c and Hell = 1. Hence 


Rc + Rc = c 


S(co)c = (I - R)c 


which on using (15) becomes 


S((i))c = -sin f))S(c)c - (1 - cos (|))s2(c)c 
But from (6) S^(c) = -I + cc^, whereas c^c = 0. Hence, 

S(oj)c = -sin (t)S(c)c + (1 - cos ()))c 
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5 . ^ (1 - cos ♦jc 

Sin (}> 

Now premultiply both sides of (Bl) by S(c) and simplify to get 

S(c)S(o))c = 2 tan - tan <t>^S(w)c 

Hence , 


c = Y S(uj)c + j cot (|)^S(c)S(w)c 

The last term in the above equation is a vector triple product, 
form of the vector triple product identity is, for any x,y, and 
S(x)S(y)z = (x^z)y - (x^y)z. Therefore 

C = S((j0)c + j cot (|)^[0J - (lo^c)c] 


The matrix 
z in E^, 
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